MODULE Commonvars

  ! GLobal MPI LANGUAGE!
  INTEGER  :: nprocs, nprocs_new, myrank, ierr, color, key, NEW_COMM, my_new_comm_rank
  EXTERNAL :: MPI_INIT, MPI_COMM_RANK, MPI_COMM_SIZE, MPI_COMM_SPLIT, MPI_FINALIZE, MPI_BARRIER, MPI_ALLREDUCE, dposv, normal_01_cdf_inv, normal_01_cdf

  ! Adult life of mayors in electoral terms
  INTEGER, PARAMETER    :: Tend = 10 !corresponds to 40 years of adult life starting from the younger mayor in the data
  INTEGER, PARAMETER    :: Retage = 8 !Period 10, they are forced to retire, not allowed to run for office
  ! Number of electoral terms simulated for each municipality for the estimation
  INTEGER, PARAMETER    :: Nsimterms = 3
  ! Number of electoral terms simulated for each municipality for the simulation of the probabilities
  INTEGER, PARAMETER    :: Nsimtermsprob = 10
  !Possible types of mayor
  INTEGER, PARAMETER    :: Ntypes = 2, Nab = 3, Nfines = 4, Nrcont = 2, Ntheta = 1, Nelap = 1
  !Probability of a particular type of mayor
  REAL(8)               :: ptype(Ntypes), pfine(Nfines), pability(Nab), pability_ud(Nab), pappeal(Nelap), pr_theta(2)
  !Nmun = number of municipalities types based on their population
  !Nactmun = number of municipalities used in the estimation
  !Ntotmun = number of municipalities used in the matrix for the nonparameteric distribution
!!$  INTEGER, PARAMETER    :: Nmun = 3, Nactmun = 2852, Ntotmunnoc = 4478, Ntotmunc = 1606
  INTEGER, PARAMETER    :: Nmun = 3, Ntotmunnoc = 4470, Ntotmunnoc_mun1 = 2364, Ntotmunnoc_mun2 = 1718, Ntotmunnoc_mun3 = 388, Ntotmunc = 1837, &
                           Nvarnonparnoctemp = 7, Nvarnonparnoc = 5, Nvarnonparc = 1, Nvarnonparsav = 3, Nvarnonparsavtemp = 4, Nactmun = 5658 !5185 !5237 !2549
  !Grids
  INTEGER, PARAMETER    :: Neduc = 2, Nsteal = 6, Naud = 2, Nwealth = 10, Nq = 3, Nfunds = 3, Npr = 1, Nmage = Retage-1, &
                           Nage = 3, Npol_q=1, Npol_age = 2, Npol_ab = 2, Npol_st = 1, MMability_steal_q = 1+Nwealth+Npol_ab+Npol_st+Npol_q, &
                           MMsteal_q = 1+Nwealth+Npol_st+Npol_q, MMage = 1+Nwealth+Npol_age, MMage_steal_q = 1+Nwealth+Npol_age+Npol_st+Npol_q
  ! Npast = 1 if audit and stealing have no effect on future wages (current version), = 3 if they have
  ! IMPORTANT FOR POLICY
  ! Npast describes the past of the current mayor: 1 = no audit; 2= audit
  INTEGER, PARAMETER    :: Npast= 2
  ! Npast_pm describes the past of a mayor that was in power more than one period ago; it is equal to 1 because their past wages do not depend on their past
  INTEGER, PARAMETER    :: Npast_pm = 1
  ! We have Npast_t = Npast+1 to take into account that voters have different wages; we use ipast_t=1 for past mayors and ipast_t=2 for voters
  INTEGER, PARAMETER    :: Npast_t = Npast_pm + 1

  !Normilizing constant
  REAL(8), PARAMETER    :: norm = 1000.0d0
!!$  !Fine functions
!!$  REAL (8)              :: fine(Nsteal) = 1.0d0*(/ (0.0d0+50000.0d0*DBLE(i) , i=0,Nsteal-1) /)/norm
  !Bound on stealing for last period
  REAL(8), PARAMETER   :: omega1 = 0d0
  REAL(8), PARAMETER   :: omega3 = 0d0
  !Number of parameters for the election probability
  INTEGER, PARAMETER    :: Nparpm = 8
  !Parameters of the election probability
  REAL(8)               :: alpha(Nparpm), p_low_ap
!!$  !Probability of Being a Mayor
!!$  REAL(8), ALLOCATABLE  :: pmayor(:,:,:,:,:,:)
  !Probability of Drawing a saving level for the first time mayor
  REAL(8)               :: pwealth(Nwealth)
  !Order of Polynomials In Value Function Approximation
  INTEGER, PARAMETER    :: Npolyst = 2, Npolysav = 2, Npolyq = 2, Npolysavc = 2
  !Size of the matrix in the polynomial used in the approximations and dVsav. The value is computed in main.f90
  INTEGER               :: Nxt
  REAL(8), PARAMETER    :: popsize(Nmun) = (/ 5000.0d0, 30000.0d0, 75000.0d0 /)
  REAL(8), PARAMETER    :: ln_popsize(Nmun) = (/ LOG(5000.0d0), LOG(30000.0d0), LOG(75000.0d0) /)
!!$  REAL(8), PARAMETER    :: popsize(Nmun) = (/ 1.0d0, 6.0d0, 15.0d0 /)

!!$  !Grid for savings
!!$  REAL(8)                :: saving(Nsav) = (/ -100000.0d0, -50000.0d0, 0.0d0, 50000.0d0, 100000.0d0, 250000.0d0, 500000.0d0, 1000000.0d0 /)/norm
!!$  REAL(8)                :: saving_sim(Nsav) = (/ 7500.0d0, 15000.0d0, 25000.0d0, 50000.0d0, 100000.0d0, 250000.0d0, 500000.0d0, 1000000.0d0 /)/norm

  ! Grid for net Wealth
!!$  REAL(8)               :: wealth(Nwealth) = (/ -500000.d0, -200000.0d0, -100000.0d0, -50000.0d0, 0.0d0, 50000.0d0, 100000.0d0, 250000.0d0, 500000.0d0, 1000000.0d0 /)/norm
!!$  REAL(8)               :: wealth(Nwealth) = (/ -500000.d0, -250000.0d0, -100000.0d0, -50000.0d0, 0.0d0, 50000.0d0, 100000.0d0, 250000.0d0, 500000.0d0, 1000000.0d0, 2500000d0, 5000000d0 /)/norm
!!$  REAL(8)               :: wealth(Nwealth) = (/ -500000.d0, -200000.0d0, -100000.0d0, -50000.0d0, 0.0d0, 100000.0d0, 250000.0d0, 500000.0d0, 1000000.0d0, 3000000.0d0 /)/norm
!!$  REAL(8)               :: wealth(Nwealth) = (/ -500000d0, -100000d0, 0d0, 50000d0, 100000d0, 250000d0, 500000d0, 1000000d0, 2500000d0, 5000000d0 /)/norm
! Grid used in estimation and policy simulations
  REAL(8)               :: wealth(Nwealth) = (/ -500000d0, -100000d0, 0d0, 50000d0, 100000d0, 250000d0, 500000d0, 750000d0, 1000000d0, 1250000d0 /)/norm
! Indeces to read in average wealth at the municipality level; iwealth_vec(1) = 4 means that average wealth for imun = 1 equals wealth(4) = 50000d0, etc
  INTEGER               :: iwealth_vec(Nmun) = (/ 4, 5, 10 /)
!!$  REAL(8)               :: wealth(Nwealth) = (/ -500000.d0, -100000.0d0, 0.0d0, 100000.0d0, 250000.0d0, 500000.0d0, 1000000.0d0, 2000000.0d0, 5000000.0d0, 10000000.0d0 /)/norm
!!$  REAL(8)               :: wealth(Nwealth) = (/ -500000.d0, -100000.0d0, 0.0d0, 100000.0d0, 500000.0d0, 1000000.0d0, 2000000.0d0, 5000000.0d0, 10000000.0d0, 20000000d0 /)/norm
!!$  REAL(8)               :: wealth(Nwealth) = (/ -500000d0, -250000d0, -100000d0, 0d0, 50000d0, 100000d0, 250000d0, 500000d0, 1000000d0, 2000000d0 /)/norm
!!$  REAL(8)               :: wealth(Nwealth) = (/ -1000000.0d0, -500000.d0, -200000.0d0, -100000.0d0, -50000.0d0, 0.0d0, 50000.0d0, 100000.0d0, 250000.0d0, 500000.0d0, 1000000.0d0, 1500000.0d0 /)/norm
!!$  REAL(8)               :: wealth(Nwealth) = (/ -500000.d0, -200000.0d0, -100000.0d0, 0.0d0, 50000.0d0, 100000.0d0, 250000.0d0, 350000.0d0, 500000.0d0, 1000000.0d0, 1500000.0d0 /)/norm
  
  !Grid for public consumption
!!$  REAL(8)               :: qconsvec(Nq) = (/ (0.0d0+50000.0d0*DBLE(i) , i=0,Nq-1) /)/norm
!!$  REAL(8)               :: qconsvec(Nq) = (/ 0.0d0, 50000.0d0, 100000.0d0, 150000.0d0, 250000.0d0, 350000.0d0, 500000.0d0 /)/norm
!!$  REAL(8)               :: qconsvec(Nq) = (/ 0.0d0, 5000000.0d0, 12500000.0d0, 25000000.0d0 /)/norm
!!$  REAL(8)               :: qconsvec(Nq) = (/ 10000.0d0, 1000000.0d0, 5000000.0d0 /)/norm
  !We use the 10th, 50th, and 90th percentile of per-capita public consumption
!!$  REAL(8)               :: qconsvec(Nq) = (/ .9269536d0, 2.395284d0, 5.840413d0 /)
!!$  !We use the 10th, 50th, and 90th percentile of log of public consumption relative to a municipality with 5000 citizens (LOG(qcons/(popsize/5000)))
!!$  !imun = 1
!!$  REAL(8)               :: qconsvec1(Nq) = (/ 8.442876d0, 9.397883d0, 10.23727d0 /)
!!$  !imun = 2
!!$  REAL(8)               :: qconsvec2(Nq) = (/ 8.349767d0, 9.132372d0, 10.49521d0 /)
!!$  !imun = 3
!!$  REAL(8)               :: qconsvec3(Nq) = (/ 8.808164d0, 9.814028d0, 10.58982d0 /)
  !We use the 10th, 50th, and 90th percentile of log of public consumption relative to a municipality with 5000 citizens (LOG(qcons/(popsize/5000)))
  !imun = 1
!!$  REAL(8)               :: qconsvec1(Nq) = (/ -.0743177d0, .8806907d0, 1.720078d0 /)
  !We compute LOG(Q/norm)=LOG(Q)-LOG(norm), where LOG(Q) is computed in the data, to have more flexibility in changing NORM
  !The three numbers correpsond to the 25th, 50th, and 75th percentiles of LOG(Q)
!!$  REAL(8)               :: qconsvec1(Nq) = (/ 14.45687d0-LOG(norm), 14.76735d0-LOG(norm), 15.08619d0-LOG(norm) /)
  !The three numbers correpsond to the 25th, 50th, and 75th percentiles of per-capital normalized public consumption
  REAL(8)               :: qconsvec1(Nq) = (/ 1.281881d0, 2.409436d0, 3.694281d0 /)
  !imun = 2
!!$  REAL(8)               :: qconsvec2(Nq) = (/ -.1674265d0, .6151783d0, 1.680982d0 /)
!!$  REAL(8)               :: qconsvec2(Nq) = (/ 15.61467d0-LOG(norm), 15.95727d0-LOG(norm), 16.34756d0-LOG(norm) /)
  REAL(8)               :: qconsvec2(Nq) = (/ 1.036066d0, 1.851396d0, 3.383907d0 /)
  !imun = 3
!!$  REAL(8)               :: qconsvec3(Nq) = (/ .2909709d0, 1.296834d0, 2.072632d0 /)
!!$  REAL(8)               :: qconsvec3(Nq) = (/ 17.11829d0-LOG(norm), 17.51419d0-LOG(norm), 18.09296d0-LOG(norm) /)
  REAL(8)               :: qconsvec3(Nq) = (/ 2.00786d0, 3.720078d0, 5.450025d0 /)
  ! If Nq = 4
  !The four numbers correpsond to the 20th, 40th, 60th, and 80th percentiles of per-capital normalized public consumption
!!$  REAL(8)               :: qconsvec1(Nq) = (/ 1.1754147d0, 1.918612d0, 2.8873839d0, 4.0755491d0 /)
!!$  REAL(8)               :: qconsvec2(Nq) = (/ .96856999d0, 1.3248601d0, 2.5028744d0, 3.8527329d0 /)
!!$  REAL(8)               :: qconsvec3(Nq) = (/ 1.7406639d0, 3.0196228d0, 4.1734753d0, 6.1039085d0  /)
  REAL(8)               :: qconsvec(Nmun,Nq)
  ! Global variable containing average qcons for mayors in the first term (potential challengers
  REAL(8), ALLOCATABLE  :: av_qcons_vec(:,:,:)
  
  ! Global variable for qconsvec for one municipality type
  REAL(8)               :: qconsvec_g(Nq)
  !Grid for public consumption
!!$  REAL(8)               :: rcont(Nrcont) = (/ (-0.5d0+0.5d0*DBLE(i) , i=0,Nrcont-1) /)
!!$  REAL(8)               :: rcont(Nrcont) = (/ 0.2704745d0, 1.307052d0, 4.96115d0 /)
!!$  REAL(8)               :: rcont(Nrcont) = (/ .5996224d0, 1.198294d0,  2.530513d0 /)
!!$  REAL(8)               :: rcont(Nrcont) = (/ (.5996224d0+1.198294d0)/2d0,  (1.198294d0+2.530513d0)/2d0 /)
  !Grid for Electoral Appeal
  !We use the name rcont because we replaced relative contributions with electoral appeal
  !fine coefficient should be 1 in the base case because we draw from a log-normal(mu,sigma)
  REAL(8), PARAMETER    :: fine_coeff = 1d0
  !Base Case: Maximum number terms allowed (2 in Brasil) or observed in the data
  INTEGER               :: Nnterms
  ! For wage policy
  REAL(8)               :: n_times_wage
  !Monthly earnings
  REAL(8)               :: wage_mayor(Nmun)
  !Monthly earnings of low-educated voters 
  REAL(8), PARAMETER    :: wage_vot(Nmun) = (/ 583.1026d0, 638.8651d0, 946.3749d0 /)/norm
  !Policy in which the mayor looses the election with probability 1 if he was caught stealing (equivalent to not allowing the mayor to run for reelection)
  INTEGER               :: flag_no_run
  !Policy in which the mayor looses the election with probability 1 if he steals more than x (equivalent to not allowing the mayor to run for reelection)
  INTEGER               :: flag_no_run_x
  !probability of being audited
  INTEGER, PARAMETER    :: Npaudits = 2
  !Audit probabilities for estimation (the ones observed in the data)
   REAL(8)              :: paudit(0:1,1:Npaudits)
  !Policy that increases the audit probability for mayors in their last term (iterm = Nnterms)
  INTEGER               :: flag_audit_nterm
  !Audit probabilities for iterm = Nnterms mayors in the flag_audit_nterm = 1 policy
  REAL(8)               :: paudit_nterm(0:1)
  !Variable used to set the audit probability in emaxfun_temp.f90, emaxfun.f90, and simulatemun.f90 before calling mayordecision.f90
  REAL(8)               :: paudit_g
  !Flag_policy_audit = 1 if we evaluate a policy in which mayors that were caught stealing are audited with probability 1
  INTEGER, PARAMETER       ::  flag_policy_audit = 0
  
  ! Index to record the results of the policy in the array where we save willingness to pay, etc
  ! Arrays needed to save the results for the base case
  REAL(8), ALLOCATABLE     :: notmayor_value_t(:,:,:), notmayor_value(:,:,:), ELNqcons_base_vec_t(:,:,:), ELNqcons_base_vec(:,:,:)
  INTEGER, PARAMETER       :: nmag = 100, nparam_subset = 4
  INTEGER                  :: param_number(nparam_subset), i_g
  REAL(8)                  :: pc_cost_prau, pr_audit_out, pc_cost_prau_out, pc_cost_wage(Nmun), n_times_wage_out, pc_cost_wage_out(3)
  REAL(8)                  :: ch_mu_a, per_ch_fine, pr_ch_age, per_ch_fine_out, prob_conv_norun, mag_ch_out(nparam_subset), prob_conv_out
  INTEGER                  :: policy_pos, ch_m_age     
  INTEGER, PARAMETER       ::              flag_policy_evaluation = 0
  INTEGER                  ::              flag_optimal_auditprob = 0
  INTEGER                  ::                    flag_base_policy = 0
  INTEGER                  ::                    flag_prau_policy = 0
  INTEGER                  ::                   flag_norun_policy = 0
  INTEGER                  ::              flag_norun_prob_policy = 0
  INTEGER                  ::                   flag_terms_policy = 0
  INTEGER                  ::                    flag_wage_policy = 0
  INTEGER                  ::                   flag_hfine_policy = 0
  INTEGER                  ::                flag_probconv_policy = 0
  INTEGER                  ::             flag_terms_norun_policy = 0
  INTEGER                  ::               flag_terms_3th_policy = 0
  INTEGER                  ::         flag_terms_norun_3th_policy = 0
  INTEGER                  ::               flag_norun_2th_policy = 0
  INTEGER                  ::        flag_terms_norun_wage_policy = 0
  INTEGER                  ::               flag_prau_wage_policy = 0
  INTEGER                  ::       flag_terms_norun_prauk_policy = 0
  INTEGER                  ::            flag_optimal_prau_policy = 0
  INTEGER                  ::       flag_term_optimal_prau_policy = 0
  INTEGER                  ::      flag_norun_optimal_prau_policy = 0
  INTEGER                  :: flag_term_norun_optimal_prau_policy = 0
  INTEGER                  :: flag_term_norun_optimal_prau_lastterm_policy = 0
  INTEGER                  ::       flag_wage_optimal_prau_policy = 0
  INTEGER                  ::       flag_optimal_wage_prau_policy = 0
  INTEGER                  ::            flag_optimal_wage_policy = 0
  INTEGER                  ::         flag_change_ability_01stdev = 0
  INTEGER                  ::  flag_change_ability_2ndterm_mayors = 0
  INTEGER                  ::          flag_change_ability_policy = 0
  INTEGER                  ::              flag_change_age_policy = 0
  INTEGER                  ::              flag_no_change_ability = 0

  !Flag_simulate_data = 1 if we simulate the data to analyze them
  INTEGER, PARAMETER       ::                  flag_simulate_data = 0

  ! The following flags determine which policy we simulate
  INTEGER                  ::                     flag_estimation = 0
 
  INTEGER                  ::                   flag_1term_policy = 0     
  !Flag_no_selec_mun = 1 if we simulate the base case without selection on municipality size  
  INTEGER                  ::                   flag_no_selec_mun = 0      
  !Flag_no_selec_funds = 1 if we simulate the base case without selection on funds
  INTEGER                  ::                 flag_no_selec_funds = 0
  !Flag_no_selec_zzpr = 1 if we simulate the base case without selection on private inputs
  INTEGER                  ::                  flag_no_selec_zzpr = 0
  !Flag_no_selec_age = 1 if we simulate the base case without selection on age
  INTEGER                  ::                   flag_no_selec_age = 0
  !Flag_no_selec_educ = 1 if we simulate the base case without selection on education
  INTEGER                  ::                  flag_no_selec_educ = 0
  !Flag_no_selec_ab = 1 if we simulate the base case without selection on ability
  INTEGER                  ::                    flag_no_selec_ab = 0
  !Flag_no_selec_sav = 1 if we simulate the base case without selection on optimal savings in the second term
  INTEGER                  ::                   flag_no_selec_sav = 0

  !flag_fixed_ability == 1 if we fix ability at some established level
  INTEGER, PARAMETER       :: flag_fixed_ability = 0
  !fixed level of ability
  REAL(8), PARAMETER       :: fixed_ability =  0.412244d0 !1.275894d0 !0.412244d0
  ! Grid for ages of mayors in Emaxfun.f90
  INTEGER               :: agegrid(Nage) = (/ 2 , 4 , 6 /)
  !Pension: income in retirement period
  REAl(8), PARAMETER    :: pension = 0.6d0
  !Bound for infeasible w.p.1
  REAL(8)               :: inf_bound_pm(Tend,Nmun,Neduc), inf_bound_pmg, inf_bound_pmgp1
  INTEGER               :: inf_bound_indexg
  !Count how many intial observations we have in which Nnterms = i
  INTEGER, ALLOCATABLE  :: countNnterms_mun1(:), countNnterms_mun2(:), countNnterms_mun3(:)
  !Matrix that contains age, education of mayors, funds, relative contributions, and savings or new mayors.
  REAL(8), ALLOCATABLE  :: nonpardistnoc_mun1(:,:,:), nonpardistnoc_mun2(:,:,:), nonpardistnoc_mun3(:,:,:), nonpardist_sav_mun1(:,:), nonpardist_sav_mun2(:,:), nonpardist_sav_mun3(:,:), nonpardistc(:), &
                           nonpardistnoc_uncond_mun1(:,:), nonpardistnoc_uncond_mun2(:,:), nonpardistnoc_uncond_mun3(:,:)
  !Nonparametric distribution for funds, private gdp, savings of new mayors, and relative contributions. Note
  REAL(8), ALLOCATABLE  :: nonpardistemax(:,:,:,:), nonpardist_fu_pr(:,:), nonpardist_wealth_fu_pr(:,:,:)
  !Nonparametric distribution for funds, private gdp, and relative contributions. Note
  REAL(8), ALLOCATABLE  :: nonpardistemaxnosav(:,:,:), nonpardistemaxnorc(:,:,:)
  REAL(8), ALLOCATABLE  :: nonpardistageduc(:,:)
  ! If we want to experiment with how funds affect stealing
  REAl(8), PARAMETER    :: per_funds = 1d0
!!$  REAl(8)               :: per_funds = 1d0
  !Grid for funds from central goverment
!!$  REAL(8)               :: funds(Nfunds) = (/ (10000.0d0+10000.0d0*DBLE(i) , i=0,Nfunds-1) /)/norm
!!$  REAL(8)               :: funds(Nfunds) = (/ 99954.44d0, 201008.8d0, 375497.8d0 /)/norm
!!$  REAL(8)                  :: funds(Nfunds) = (/ 58864.41d0, 120766.2d0, 234071.6d0 /)/norm
!!$  REAL(8)               :: funds(Nfunds) = per_funds*(/ 978340.8d0, 1326480.0d0, 2070997.0d0 /)/norm
  REAL(8)               :: funds(Nfunds) = per_funds*(/ 669618.8d0, 1018890d0, 1803933d0 /)/norm
  !Grid for stealing
  REAL(8)               :: fr_stoleni(Nsteal) = (/ 0.0d0, 0.025d0, 0.05d0, 0.1d0, 0.2d0, 0.5d0 /)
  REAL(8)               :: stealing(Nsteal)

  REAL(8),PARAMETER     :: steal_bound1 = 1d0/norm, steal_bound2 = 75000d0/norm, steal_bound3 = 11111125000d0/norm

  !Grid for private inputs in the production function for public consumption
!!$  REAL(8)               :: zzpr(Npr) = (/ (10000.0d0+10000.0d0*DBLE(i) , i=0,Npr-1) /)/norm
!!$  REAL(8)               :: zzpr(Npr) = (/ 8874.912d0, 18129.54d0, 42117.92d0 /)/norm
!!$  REAL(8)               :: zzpr(Npr) = 1000.0d0*(/ 4485.7d0, 9645.2d0, 22747.3d0 /)/(norm*1000d0)
!!$  REAL(8)               :: zzpr(Npr) = (/ 8982518.0d0, 18400000.0d0, 44600000.0d0 /)/norm
!!$  REAL(8)               :: zzpr(Npr) = (/ 2.129d0, 2.587d0, 3.132d0 /)
  REAL(8)               :: zzpr(Npr) = (/ 2.553d0 /) ! zzpr(Npr) = (/ 2.142d0, 2.553d0, 3.133d0 /)
  !Audit = 1 if the municipality was audited, 0 otherwise
  REAL(8)                  :: audit(Naud) = (/ 0.0d0, 1.0d0/)
  !Annual Labor supply if never a mayor in months
  REAL(8)                  :: lsnotmayor = 12.0d0
  !Annual Labor supply if mayor in months
  REAL(8)                  :: lsmayor = 12.0d0
  !Annual Labor supply if mayor in the past in months
  REAL(8)                  :: lspastmayor = 12.0d0
  !Risk-free interest rate
  REAL(8), PARAMETER       :: R = 1.05d0
  !discount factor
  REAL(8), PARAMETER       :: betadf = 0.98d0
  !Number of parameters in the production function for public consumption
  INTEGER, PARAMETER       :: Nparq = 9
  !parameters of the production functions for public consumption
  REAL(8)                  :: alphaq(Nparq)
  ! Relative preferences for public consumption parameters
  REAL(8)                  :: theta
  ! Electoral appeal parameters
  REAL(8)                  :: alpha_ie, alpha_ea, zeta, alpha_st, alpha_q, alpha_nst
  !population polynomial for production function
  REAL(8)                  :: polypop(Nmun)
  !Number of parameters in the wage function for past mayors
  INTEGER, PARAMETER       :: Nwage = 8
  !parameters of the wage function for past mayors
  REAL(8)                  :: alphaw(Nwage), sigmaw
  !gamma and delta are the preference parameters for private and public consumption
  REAL(8)                  :: gamma, delta
  !eta determines the degree of rivalry of public consumption in the utility functino
  REAL(8)                  :: eta
  !Extra Utility agents get from being in power
  REAL(8)                  :: upower
  !Cost of running parameters: fraction and level
  REAL(8)                  :: cost_run_frac
  !Standard deviation of the private input shocks and stealing, and parameters that reduces the fine for some types of mayors
  REAL(8)                  :: sigma_zzpr, sigma_run, mu_steal, sigma_steal, frac_fine_g
  !Parameters of the continuous distribution of fines
  REAL(8)                  :: mu_fine, sigma_fine, delta_fine
  !Vector characterizing the discretized fraction of the fine, if Nfines = 3, then (/ 1/4, 2/4, 3/4 /)
  REAL(8)                  :: frac_fine(Nfines)
  !Probability of being convicted
  REAL(8)                  :: prob_conv
  !Parameters of the continuous distribution of fines
  REAL(8)                  :: mu_a, sigma_a
  !Vector characterizing the discretized value of ability
  REAL(8)                  :: a_discr(Nab)
  !grip point for savings where savings = 0.0
  INTEGER                  :: w0
  !global variables for mayor and not mayor decisions
  REAL(8)                  :: waget, fundst, zzprt, savingt, savingct, ability_g, finet, stealt,  ELNqqt, agemayort, totrest, boundt, polyt, popsizeg, ipaug, &
                              min_ability, av_saving, av_zzpr, av_funds, thetai_g, newsteal_g, new_fr_stolen_g, past_g
  REAL(8)                  :: vemaxoptg, pconsoptg, qconsoptg, newsavoptg, stealoptg, utiloptg, emaxoptg, zzpuoptg, vemaxoptpmg, pconsoptpmg, newsavoptpmg, utiloptpmg, emaxoptpmg, &
                              notmayor_low_value_g, pmwage_g, med_wealth_g, ELNqcons_high_g, ELNqcons_g
  INTEGER                  :: imunt, iedut, ipastg, itermg, MMg, d_past_emax, timeg, med_age_g, itermi_g
  !number to assign to util if -infinity
  REAL(8), PARAMETER       :: bignumber = -10000000000.0d0, check_bignumber = -1000000000.0d0, focbignumber = 10000000000.0d0, check_focbignumber = 1000000000.0d0 
  !parameters of the value function
  REAL(8), ALLOCATABLE     :: betat(:),betat1(:),betat2(:),betat3(:)
  !parameters used in the approximation of the value functions
  REAL(8), ALLOCATABLE     :: betanewcr(:,:,:,:,:,:,:,:), betapmayor(:,:,:,:,:,:,:,:,:,:,:), betapmayortemp1(:,:,:,:,:,:,:,:,:,:,:,:), &
                              betapmayornor(:,:,:,:,:,:,:), betasimrun(:,:,:,:,:,:,:,:), betanornew(:,:,:,:,:,:,:), &
                              beta_voter_i(:,:,:,:,:,:,:,:,:), beta_voter_c(:,:,:,:,:,:,:)






  !Number of simulations for each municipality
  INTEGER, PARAMETER       :: Nsim = 100 ! WITH NSIM = 500 THE RESULTS ARE ABOUT THE SAME







  !Number of variables in simulation, in the simulation for the initial probabilities, and in the matrix of initial conditions
  INTEGER, PARAMETER       :: Nvar = 36, Nvarprob = 4, Nvarinit = 18
  INTEGER                  :: Nvarreg, Ninst
  !Data for the initial conditions
  REAL(8), ALLOCATABLE     :: initdata(:,:)
  !Variable to stop program
  INTEGER, PARAMETER       :: istop = 0
  !Number of parameters to be estimated
  INTEGER, PARAMETER       :: nparam = 11
  !Factors multiplying the parameters to make all of them between 0 and 1
  REAL(8)                  :: factor_param(nparam)
  !Number of moments required to construct the moments used in the estimation
  INTEGER, PARAMETER       :: Nmomtemp = 17
  !Number of moments used in the estimation
  INTEGER, PARAMETER       :: Nmom = 11
  !Number of moments in data
  INTEGER, PARAMETER       :: Nmomdata = 11
  !Vector of parameters used in standarderrors_SMM to find the parameters that minimize fval
  REAL(8)                  :: opt_param(nparam), dgdp_g(Nmom,nmag)
  !flag used in main and standarderrors_SMM to find the parameters that minimize fval
  INTEGER                  :: flag_lowerfval
  !Vector of simulated and data moments used in the estimation
  REAL(8)                  :: simdatamom(Nmom,1), simmom(Nmom,1), datamom(Nmomdata,1), transform(Nmom,Nmomdata), norm_factor(Nmomdata)
  !Vector of data moments and its covariance matrix
  REAL(8)                  :: data(Nmomdata,1), dataV(Nmomdata,Nmomdata)
  !Variance-covariance matrix and its inverse
  REAL(8)                  :: V(Nmom,Nmom), INV_V(Nmom,Nmom)
  !Total number of parameters
  INTEGER, PARAMETER       :: nparamsim = 25 + 2*Nmun
  !Bounds for parameters
  REAL(8), ALLOCATABLE     :: ubp(:), lbp(:)
  !If iprint = 1 then print
  INTEGER, PARAMETER       :: iprint = 1
  !Flag_sim = 1 if we simulate the data to compute the data moments
  INTEGER, PARAMETER       :: flag_sim = 0
  !flag_sim_data = 1 inside moments.f90 when we simulate the data 
  INTEGER                  :: flag_sim_data = 0
  !Random draws to generate simulated data
  REAL(8), ALLOCATABLE     :: udraw(:,:,:,:), utypeprob(:,:), fine_norm(:,:,:), ability_norm(:,:,:), theta_norm(:,:,:), age_udraw(:,:,:,:)
  !Draws from a random normal to generate simulated data
  REAL(8), ALLOCATABLE     :: rnorm(:,:,:), rnormprob(:,:,:), zzpr_shock(:,:,:)
  INTEGER                  :: flag_here = 0, flag_simulation = 0, flag_Tend = 0, flag_fafb
  !Subsistence level of utility in case there is no public or private consumption
  REAL(8), ALLOCATABLE     :: subutil(:)
  !Number close to zero
  REAL(8), PARAMETER       :: small_no = 0.01/norm  !1.0d-5
  !Subsistence level of public or private consumption in case there is no positive public or private consumption
  REAL(8), PARAMETER       :: subpcons = 1000d0/norm
  REAL(8), PARAMETER       :: subqcons_vec(Nmun) = 10000d0*(/ 1d0/(norm*popsize(1)), 1d0/(norm*popsize(2)), 1d0/(norm*popsize(3)) /)
  REAL(8)                  :: subqcons
  !Starting value for lower bound in zeroin
  REAL(8), PARAMETER       :: start_a = 0.0d0
  !Slope parameter for dVst at high levels of savings
  REAL(8), PARAMETER       :: slopedVst = 1.0d-3
!!$  !Lower bound to be used in maximinations
  REAL(8), PARAMETER       :: lbound = -10000000.0d0/norm
  REAL(8), PARAMETER       :: lbound_sav = -1d+8
  !Lower bound to compute the initial conditions in maximinations
  REAL(8), PARAMETER       :: lboundit = -90000000000.0d0/norm, infty = 1d+8, check_infty = 1d+7, min_concave = -1d-6, min_positive = 1d-4, min_negative=-1d-4, min_qcons = 1d-6, inc_slope = 0.25d0, &
                              per_min_util = 0d0, per_min_util_steal = 0d0, ini_xstar_steal = 100001d0/norm, ini_xstar_sav = 50001d0/norm, sub_level = 0d0/norm, min_slope = -1d-6, Rsquared_lb = 0.5d0, &
                              inc_wealth = 1d0, l_steal_prob = 0d0, ch_slope = 1d0, max_frac_fine = 10d0
  REAL(8)                  :: l_fr_stolen = 0d0, l_fr_stolen_x = 0d0
  REAL(8), ALLOCATABLE     :: probemax1(:,:,:,:,:,:,:,:,:,:,:,:,:,:),probemax2(:,:,:,:,:,:,:,:,:,:,:),probemax3(:,:,:,:,:,:,:,:,:,:,:,:),probemax4(:,:,:,:,:,:,:,:,:,:,:,:),probqcons_c(:,:,:,:,:,:,:,:), probqcons_i(:,:,:,:,:), Pr_i_win(:,:,:,:,:,:,:,:,:,:,:,:,:)
  REAL(8), ALLOCATABLE     :: yp(:), xp(:,:)
  REAL(8)                  :: accrcy, machine_eps, previous_fvec(2), previous_xx(2), min_utilg(Neduc,Nmun), av_util_mayors_g(Nmun,Tend), &
                              ub_pr, lb_pr, ub_gamma, lb_gamma, lb_a1, ub_a1, lb_a2, ub_a2, lb_a3, ub_a3, lb_a7, ub_a7
  INTEGER                  :: size_probit, previous_info
  ! number of iterations
  INTEGER                  :: iterations
  ! Minimization method: if flag_sa = 1 we use simulated annealing, if flag_sa = 0 we use Nelder-Mead                          
  INTEGER                  :: flag_sa = 0
  ! Minimization method: if flag_sopt_sa = 1 we use simulated annealing, if flag_sopt_sa = 0 we use Nelder-Mead
  INTEGER :: flag_sopt_sa
  INTEGER :: flag_nan_minim
  INTEGER, PARAMETER       :: flag_drawrandom = 0
  INTEGER, PARAMETER       :: flag_drawparam = 0
  ! flag_eff_weight = 1 if we use the efficient weighting matrix
  INTEGER, PARAMETER       :: flag_eff_weight = 1
  INTEGER, PARAMETER       :: flag_compstderrs = 0
  ! flag_bootstrap = 1 if we estimate the parameters using one of the boostrap samples
  INTEGER, PARAMETER       :: flag_bootstrap = 0
  ! nboot = one of bootstrapped estimations (1-200)
  INTEGER                  :: nboot, first_proc
  ! Nboot_min = smallest bootstrap in loop, Nboot_max =highest bootstrap in loop; nprocs divided by Nboot_max-Nboot_min+1 should give 8 if we want 8 procs for each bootstrap (64/8=8, 128/8=16)
  INTEGER, PARAMETER       :: Nboot_min = 1, Nboot_max = 1
  REAL(8), ALLOCATABLE     :: bootstrap_vect(:,:,:), bootstrap_vec(:,:,:)
  ! cost_policy = 1 when we account for the cost of audits or increased wages
  INTEGER                  :: cost_policy
  ! Number of production function parameters estimated outside the model: Nparq + 2 (mean and std deviation of abilty)
  REAL(8)                  :: prod_func_param(Nparq+2)
  !flag_trunc_norm = 1 if we use a truncated normal to model the fine distribution, = 0 if we use a beta distribution
  INTEGER, PARAMETER :: flag_trunc_norm = 1
  INTEGER, PARAMETER :: add_draws = 200000
  ! To be used in probit.f90
  REAL(8),PARAMETER     :: sqr2pi = 1d0/SQRT(3.1415926535897932384626*2d0)
  ! Probability of observing a type of education in the data: edu = 1 (high school or less); = 2 (more than high school)
  REAL(8)               :: probedu(Neduc) = (/ 0.6192d0, 0.3808d0 /)

!qrsh -l i,h_rt=4:00:00
!qrsh -l h_rt=11:00:00,h_data=8G
!module load stata; cd /u/home/m/mmazzoc/project/MyPapers/Mayors/Estimation/Stata; xstata-se &

END MODULE Commonvars
